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Abstract 


Miscanthus lutarioriparius is an endemic species that grows along the middle and lower reaches of the Yangtze 
River and is a valuable source of germplasm for the development of second-generation energy crops. The plant 
that propagates via seeds, stem nodes, and rhizomes shows high phenotypic variation and strong local adapta- 
tion. Here, we examined the magnitude and spatial distribution of genetic variation in M. lutarioriparius across 
its entire distributional range and tested underlying factors that shaped its genetic variation. Population genetic 
analyses were conducted on 644 individuals from 25 populations using 16 microsatellite markers. M. lutarioripa- 
rius exhibited a high level of genetic variation (Hg = 0.682-0.786; A, = 4.74-8.06) and a low differentiation 
(Fsr = 0.063; Dest = 0.153). Of the total genetic variation, 10% was attributed to the differences among popula- 
tions (df = 24, P < 0.0001), whereas 90% was attributed to the differences among individuals (df = 619, 
P < 0.0001). Genetic diversity did not differ significantly across longitudes and did not increase in the popula- 
tions growing downstream of the Yangtze River. However, significant associations were found between genetic 
differentiation and spatial distance. Six genetic discontinuities were identified, which mostly distributed among 
downstream populations. We conclude that anthropogenic factors and landscape features both contributed to 
shaping the pattern of gene flow in M. lutarioriparius, including long-distance bidirectional dispersal. Our results 
explain the genetic basis of the high degree of adaptability in M. [utarioriparius and identify potential sources of 
new germplasm for the domestication of this potential second-generation energy crop. 
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Introduction 


Analysis of population connectivity plays an important 
role in protecting endangered species and managing 
economically important species, which was defined as 
the degree to gene flow affecting evolutionary processes 
within populations (Hedgecock et al., 2007; Lowe & Al- 
lendorf, 2010). Likewise, plants’ genetic diversity can 
contribute significantly to their adaptability to environ- 
mental changes such as global warming, drought, and 
urbanization (Hughes et al., 2008; Medrano & Herrera, 
2008; Markert et al., 2010; Verhoeven et al., 2011; Wang, 
2012; Raymond et al., 2013). Therefore, understanding 


Correspondence: Jianqiang Li, tel. +86 27 87510330, fax: +86 27 
87510251, e-mail: liiq@wbgcas.cn; Tao Sang, tel. +86 10 62836172, 
fax: +86 10 62590843, e-mail: sang@ibcas.ac.cn 


gene-flow patterns and genetic diversity could facilitate 
conservation practices and exploitation of genetic 
resources for utilization. 

Traditionally, methods like Wright's FST-statistic 
(1969) and Slatkin’s private allele have been used to 
analyze gene flow indirectly by measuring the spatial 
distribution of gene frequencies over time. However, 
these methods are inadequate to capture short-term 
changes in movement patterns, including populations 
that are not in equilibrium (Abdo et al., 2004; Paetkau 
et al., 2004). With advanced molecular tools, a landscape 
genetic approach was developed to study the interac- 
tion between landscape features and gene flow, which 
enables characterizing population genetic processes of 
plants inhabiting river systems, for instance (Manel 
et al., 2003; Pollux et al., 2009; Murphy et al., 2010; 
Andrew et al., 2012; Manel & Holderegger, 2013; Rico 


© 2015 The Authors. Global Change Biology Bioenergy published by John Wiley & Sons Ltd. 
This is an open access article under the terms of the Creative Commons Attribution License, 
which permits use, distribution and reproduction in any medium, provided the original work is properly cited. 1 


2 J. YAN et al. 


et al., 2014). This approach can be used to explain the 
spatial genetic variation patterns in a species by charac- 
tering the clines, isolation by distance, and genetic 
boundaries to gene flow among populations and by 
relating theses genetic structures to landscape features 
(Meeuwig et al., 2010; Song et al., 2013; Wei et al., 2013; 
Taylor & Hoffman, 2014). 

River systems provide an opportunity to study the 
interactions between genetic variation and eco-geo- 
graphical factors of hydrochloric (i.e., water dispersing) 
plants due to the unidirectional flow of water (Bergl & 
Vigilant, 2007; Honnay ef al., 2010; Orsini et al., 2013; 
Catford et al., 2014). Generally, hydrological dispersion 
results in the increased genetic diversity along the river 
systems by seeds and clones from flood damaged 
stands-washing the steams and rhizomes downstream. 
This hypothesis was confirmed for Myricaria laxiflora in 
the Three Gorges Valley of the Yangtze River (Liu et al., 
2006), Sparganium emersum in a stretch of Nier River 
(Pollux et al., 2009) and Heliconia metallica along a 
stretch of river in the Amazon lowlands (Schleuning 
et al., 2011), but other studies have not found the 
expected pattern (Chen et al., 2009; Wei et al., 2013). 
Thus, understanding the patterns of gene flow and pop- 
ulation connectivity of the riparian plant is valuable for 
studying the genetic variation and adaptive potential of 
the species. 

Miscanthus lutarioriparius is a perennial grass that pri- 
marily grows in wetlands along the middle and lower 
reaches of the Yangtze River (29-31°N; 112-119°E). It 
was recognized as a separate species from Miscanthus 
sacchariflorus according to the thicker and taller tillers 
and vigorous growth along seasonally flooded river 
banks and lake shores in central China (Chen & Renvo- 
ize, 2006), despite some remaining taxonomic contro- 
versy (Sun et al., 2010). Miscanthus lutarioriparius can 
propagate clonally by rhizomes and sexually through 
outcrossing. The hairy caryopsis seeds are capable of 
long-distance dispersal through biotic and abiotic dis- 
persal agents. Because of its commercial utilization by 
the paper industry, the propagules (seeds and rhi- 
zomes) are migrated long distances by humans through 
water and land transportation; especially, Miscanthus 
species are considered to be promising second-genera- 
tion energy crops (Jørgensen & Schwarz, 2000; Clifton- 
Brown et al., 2004; Somerville et al., 2010; Sang & Zhu, 
2011; Clark et al., 2014). Some recent studies of Miscan- 
thus genetic diversity were carried out, including 
Miscanthus x giganteus, M. sacchariflorus, and Miscanthus 
sinensis (Chae et al., 2014; Clark et al., 2014, 2015; Głowa- 
cka etal., 2015). To date, little is known about the 
genetic diversity of M. lutarioriparius across its distribu- 
tion range, although it is particularly promising energy 
crop with high biomass production and water-use effi- 


ciency, and strong ability to grow in the stressful envi- 
ronments (Sang, 2011; Liu et al., 2012, 2014; Yan et al., 
2012, 2015; Liu & Sang, 2013; Mi et al., 2014; Xu et al., 
2015). 

We used 16 highly polymorphic microsatellite mark- 
ers to analyze 25 M. Iutarioriparius populations across its 
native distribution range. The aims of our study were to 
examine the magnitude and spatial distribution of 
genetic variation and identify the gene-flow patterns in 
M. lutarioriparius. For these purposes, we first estimated 
the genetic diversity and population connectivity of the 
surveyed populations. Based on landscape features and 
eco-geographical patterns, we predicted that M. lutarior- 
iparius is subject to substantial genetic recombination 
and bidirectional migration. Our study sheds light on 
the genetic structure and gene flow of a riparian plant 
and identifies populations which may be particularly 
valuable sources of new germplasm for energy crop 
development. 


Materials and methods 


Sample collection 


From 2008 to 2010, a total of 644 of M. lutarioriparius individu- 
als were collected from 25 populations in the middle and lower 
reaches of Yangtze River (Fig. 1), which covered the entire 
native range of the species. Nineteen populations were col- 
lected along the Yangtze River, and six populations were col- 
lected in land. From each plant, cleaned plant leaves 
(2 x 5 cm) were sampled every 10-20 m across linear transects 
and immediately stored in plastic bags filled with silica gel for 
subsequent DNA extraction. 


Microsatellite genotyping 


Plant genomic DNA was isolated from approx. 1 g of dried 
leaves using a modified cetyltrimethy lammonium bromide 
(CTAB) protocol (Doyle & Doyle, 1987). DNA quality was 
tested on 0.8% agarose gels. After measuring DNA concentra- 
tion with an Eppendorf BioPhotometer, samples were diluted 
to 5 ng ult. Twenty-five individuals (one per population) 
were chosen for microsatellite primer selection. A total of 35 
microsatellite loci (Hung et al., 2009; Zhou et al., 2011, 2013) 
were tested in the study, and sixteen microsatellite loci were 
amplified by PCR for further analysis (Table S1, Supporting 
information). PCRs were carried out in a GeneAmp PCR Sys- 
tem 9700 (Applied Biosystems, Foster City, CA, USA) with 
10 uL volumes containing 1* PCR buffer, 1.5 mm MgCl2, 
200 um each dNTP, 0.4 uM fluorescently labeled forward primer 
(FAM and HEX), 0.4 um reverse primer, and 0.5 U Expand 
High-Fidelity Polymerase. PCR conditions were as follows: 
94 °C for 5 min, 32 cycles of 94 °C for 50 s, 45 s of annealing at 
54-60 °C, 1 min of elongation at 72 °C, a final extension step of 
72 °C for 10 min. Microsatellite products were run on an ABI 
3730XL automated sequencer at the Sangon Biotech and geno- 
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Fig. 1 Locations of the sampled Miscanthus lutarioriparius populations across its distribution range. Population identities correspond 
to codes in Table 1. Red dots indicate samples along the Yangtze River, and green dots indicate samples from land in the middle and 


lower Yangtze River. 


typed using GeneMarker 1.95. Based on the results of flow 
cytometry (Li et al., 2013) and the peaks of reads, we elimi- 
nated the probable tetraploid individuals. A total of 644 indi- 
viduals were included in the dataset for analysis (ranged from 
12 to 33 individuals and mean 25.8 individuals/population). 
Accessions were scored at each locus for two alleles, corre- 
sponding to the two largest peaks shown by GeneMarker 1.95. 
Potential scoring errors and null alleles were checked by micro- 
CHECKER v2.2.3 (Van Oosterhout et al., 2004). 


Data analysis 


Genetic diversity analysis 


We used GenAIEx 6.0 (Peakall & Smouse, 2006) to iden- 
tify the clonal structure of each population. Repeating 
genotypes within populations were excluded from all 
further analyses. Next, we used GENEPOP v4.0.11 (Ray- 
mond & Rousset, 1995) to test whether loci and popula- 
tions were in Hardy—Weinberg equilibrium (HWE), 
using Fisher’s exact test (Guo & Thompson, 1992), and 
ran a Markov chain for 10 000 steps following 10 000 
dememorization steps for a standard Bonferroni correc- 
tion. We also used FSTAT version 2.9.3 (Goudet, 1995) 
to analyze linkage disequilibria (LD) for all pairs of loci 
in each population, the number of alleles (A), observed 
heterozygosity (Ho), genetic diversity (Hg), fixation 
index (Fıs) averaged across loci for each populations, as 


well as genetic differentiation among populations (Fsr) 
and allelic richness (A,). To test the unidirectional 
genetic diversity hypothesis, observed heterozygosity 
(Ho) and allele richness (A,) of the populations were 
correlated with longitude of the population origin using 
a Spearman rank correction with spss 16.0 (SPSS Inc., 
Chicago, IL, USA). The Yangtze River is generally 
divided into three sections: the upper reaches, from the 
river’s source to Yichang (more than 4500 km long), the 
middle reaches, from Yichang to Hukou (29.71°N, 
115.97°E), (950 km long), and the lower reaches, from 
Hukou to the river mouth at about 930 km long (Yang 
et al., 2006). Based on the geographic boundaries of the 
Yangtze River, we divided populations as upstream and 
downstream with the boundary point at longitude 
115.7°E and calculated the correlation between genetic 
variation (Ho and A,) and longitude using spss 16.0. 


Genetic differentiation and population structure 


Genetic structure of the species was investigated using 
a number of approaches: (i) analysis of molecular 
variance (AMOvA) and F-statistics, (ii) PCoA, and (iii) 
nonhierarchical Bayesian clustering. First, AMOVA was 
performed on all individuals using GenAlEx 6.0 
(Peakall & Smouse, 2006). Pairwise F-statistics (Weir & 
Cockerham, 1984) were calculated to estimate the 


© 2015 The Authors. Global Change Biology Bioenergy published by John Wiley & Sons Ltd., GCB Bioenergy, doi: 10.1111/gcbb.12278 


4 J. YAN et al. 


genetic differentiation between populations, based on 
the differences in allele frequencies. Their significance 
levels were estimated with 6000 permutations as 
implemented in rstaT 2.9.3 (Goudet, 1995). Because Fer 
is dependent on marker variability (Hedrick, 2005; 
Jost, 2008), we additionally calculated pairwise Dest 
(Jost, 2008), an alternative measure of genetic differen- 
tiation, using the program SMOGD (Crawford, 2010). 
Second, we performed principal coordinate analysis 
(PCoA) at the different levels based on Nei’s genetic 
distance (Na, Nei et al., 1983). Third, to explore the 
genetic groups within the samples, nonhierarchical 
Bayesian clustering was performed using the program 
STRUCTURE 2.3.4 (Pritchard et al., 2000). The program 
was used with a burn-in of 100 000 replicates, fol- 
lowed by 1000 000 replicates of the Markov chain 
Monte Carlo (MCMC) simulation, with the admixture 
model and the assumption of correlated allele frequen- 
cies among populations. Three iterations were per- 
formed for each value of K from 1 to 15. The output 
was interpreted with Structure Harvester (Earl & 
vonHoldt, 2012) using the methods of Pritchard et al. 
(2000) and Evanno etal. (2005) and visualized by 
DISTRUCT v. 1.1 (Rosenberg, 2004). 


Population genetic structure landscape 


Mantel tests were performed in IBD 2.0 (Bohonak, 
2002) with 10000 permutations between pairwise 
genetic distance ([Fsr/(1—Fsr)]; [Dest/(1—Dest)]) and 
geographic distance (km) in all populations. Addition- 
ally, spatial autocorrelation analysis was performed 
on genetic distances between all individuals over 
multiple populations using GenAIEx 6.0 (Smouse & 
Peakall, 1999; Peakall & Smouse, 2006). Multiple dis- 
tance classes were tested to have the best chance of 
identifying spatial genetic structure, and these analy- 
ses were repeated at variable distance classes and 
even sample classes. Permutation and bootstrap tests 
were conducted with 10 000 replicates in each case. 
Barrier 2.2 (Manni et al., 2004) implementing Monmo- 
nier’s algorithm was used to identify the geographical 
boundaries associated with genetic differentiation 
among populations. Barriers calculated geographic 
discontinuities in the genetic structure according to a 
multilocus genetic differentiation (Fsr; Weir & Cocker- 
ham, 1984), along with their geographical coordinates 
for each population. Analyses were conducted both 
for each of the sixteen loci separately, and also for 
the complete dataset. Only barriers supported by at 
least half of the loci (i.e., >8 loci) were allowed as 
‘consensus barriers’. The analysis of individual loci 
estimated how many loci supported each identified 
barrier. 


Gene-flow and population history 


We used GENECLASS 2.0 (Piry et al., 2004) to perform the 
assignment test with the Bayesian method (Rannala & 
Mountain, 1997) and to detect the first-generation 
migrants among populations with the frequency-based 
method (Paetkau et al., 1995). The assignment test was 
used to assign or exclude the possible origins of individ- 
uals in reference populations, with 1 000 000 simula- 
tions and an alpha level of 0.01. The first-generation 
migrants among populations were conducted with 
Monte Carlo resampling of 1000 individuals and a 
threshold of 0.01, which identified the immigrant indi- 
viduals that were not born in the population. The migra- 
tion rate could indicate the gene flow at a particular 
time, but it may not accurately reflect long-term patterns 
of gene flow (Slatkin, 1987). Fsr and private allele meth- 
ods were used to estimate historical rates of gene flow 
(Nm). Fsr method was executed according to Wright’s 
island model of population genetic structure (Fsr~1/ 
(1 + 4Nm)) (Slatkin & Barton, 1989). The private allele 
method was used to estimate the effective numbers of 
migrants per generation using GENEPOP v4.0.11 (Ray- 
mond & Rousset, 1995). BOTTLENECK 1.2.02 (Piry et al., 
1999) was used to assess whether a population exhibited 
a significant bottleneck. The program assumes that the 
population after a bottleneck shows a larger heterozy- 
gosity than what is expected under the mutation-drift 
equilibrium. We conducted the Wilcoxon sign-rank test 
with a two-phase mutational model (TPM, Di Rienzo 
et al., 1994) and two other mutational models (stepwise 
mutation model and infinite allele model). The TPM 
was set at 90% single-step mutations, assuming a con- 
servative variance among multiple steps of 10. With the 
expectation that rapid reduction in allelic diversity 
under a bottleneck would cause a mode-shift distortion 
(Luikart et al., 1998), the allele frequency distribution 
test was also performed to examine frequencies of all 
alleles in populations and compare them with the distri- 
bution expected at mutation-drift equilibrium. 


Results 


Defining genotypes and estimating genetic diversity of 
Miscanthus lutarioriparius 


Among the multilocus genotypes, three pairs of clones 
were identified: plant #3 was similar to plant #4, #28 was 
similar to plant #29 (both within LU8), and plant #1 was 
similar to #2 in LU16. We calculated the maximal clone 
sizes of the species with 13.153 m and 21.213 m in LU8 
and LU16, respectively. The subsequent analyses were 
performed on the 641 genetically distinct individuals. As 
allele frequencies per locus per population were not 
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consistent in their pattern for other loci or populations, 
all sixteen loci were used for subsequent analyses. A 
total of 203 alleles were identified, ranging from 5 (C-1) 
to 19 (F-1), with an average of 12.68 alleles per locus 
(Table S1). The mean number of alleles within popula- 
tions ranged from 5.125 (LU8) to 8.875 (LU10) (Table 1). 
The mean observed heterozygosity (Ho) and expected 
heterozygosity (Hg) per population were 0.866 and 
0.766, respectively (Table 1). The smallest and the high- 
est Hy values were observed in the LU8 and LU13 popu- 
lations, which were growing along the tributary of the 
Yangtze River. After excluding populations with less 
than 20 individuals, allelic richness ranged from 4.74 
(LU8) to 8.06 (LU10) (Table 1). The average value of the 
fixation index was very low (mean: —0.161, range: —0.36 
to —0.07; Table 1), suggesting that M. I[utarioriparius is 
strictly outcrossing, and our samples were derived from 
randomly mating populations (panmixis). 

Significant departures (P < 0.05) from Hardy—Wein- 
berg equilibriums were detected in 170 of the 400 
population-loci comparisons. According to all loci and 
populations, this result may arise from the low fre- 
quency of null alleles or small sampling effects. 


Genotypic disequilibrium was tested in 120 pairs of loci 
combinations by Fisher’s exact test. No pairs of loci 
showed a significant disequilibrium after Bonferroni 
correction, which revealed that loci were physically 
unlinked and statistically independent. 

Spatial variation in genetic diversity was examined 
across the species’ range. No significant spatial patterns 
in Ho and A, among populations at longitude were 
detected (Ho: r= —0.005, P = 0.998; A,: r= —0.140, 
P = 0.544; Fig. 2). Upstream populations (longitude 
breakpoint at 115.7°E) showed a comparable average Hg 
(0.744) and A, (7.323) as downstream population (0.748 
and 6.996, respectively). Linear correlations revealed a 
nonsignificant increase for Ho (r= 0.160, P = 0.553; 
Fig. 2a) and A, (r = 0.192, P = 0.476; Fig. 2b) in upstream 
population. In downstream populations, genetic diver- 
sity showed a trend for positive correlation with Ho 
(Fig. 2a) and for negative correlation with A, (Fig. 2b). 


Genetic differentiation and population structure 


Table S2 shows a matrix of the pairwise population 
differentiation estimated by Fsr (below the diagonal) 


Table 1 Summary measures of genetic diversity for each population, including the location of populations, sample size (N), average 
of allele numbers (A), observed heterozygosity (Ho), expected heterozygosity (Hg), allelic richness (A,), and inbreeding coefficient 
(Fis). Bold populations are from the land and others are along Yangtze River 


Populations Latitude (N) Longitude (E) N A Ho Hg A, Fis 
LU1 30°24'09.6" 111°41'46.7" 31 8.125 0.843 0.743 7.42 —0.134 
LU2 28°53'29.4" 112°01/32:0” 29 8 0.842 0.706 73 —0.193 
LU3 28°53'/24.5" 112°06'49.2” 30 7.6875 0.879 0.741 6.93 —0.185 
LU4 28°52'13.2" 112°22'26.3" 27 7.625 0.824 0.752 7.16 —0.095 
LU5 29°46'27.7" 112°36'11.7" 30 8.5 0.902 0.744 7.54 —0.212 
LU6 29°09'51.8" 112°44’06.3”" 29 8.5 0.831 0.725 7.59 —0.146 
LU7 29°47'51.4" 112°51'40.1" 29 8.5 0.823 0.744 7.81 —0.105 
LU8 30°27'04.8” 112°56'02.3” 33 5.125 0.927 0.682 4.74 —0.36 
LU9 29°24'48.8" 113°04’05.1" 30 8.375 0.903 0.748 7.48 —0.207 
LU10 29°52'02.1" 113°34’44.8" 30 8.875 0.841 0.752 8.06 —0.118 
LU11 29°57'36.6" 113°40'48.8”" 30 7.4375 0.824 0.726 6.86 —0.133 
LU12 30°09/13.3” 113°42'32.3" 24 8.0625 0.907 0.777 7.71 —0.17 
LU13 30°39'48.9" 113°43/35.7" 30 8.6875 0.921 0.786 7.9 —0.17 
LU14 30°08'24.9" 114°04'18.3" 30 8.5 0.925 0.774 7.69 —0.194 
LU15 30°23/56.9" 115°00'26.6” 23 8.0625 0.882 0.754 7.75 —0.169 
LU16 30°13'06.8" 115°01'08.3" 30 7.75 0.817 0.764 7.23 —0.07 
LU17 29°44'38.65" 116°09'55.49" 12 6.437 0.867 0.756 — —0.15 
LU18 29°59'25.88" 116°39'08.05" 28 7.687 0.825 0.726 7 —0.135 
LU19 30°29'38.40" 117°00'07.34" 29 8.25 0.898 0.7665 7.51 —0.171 
LU20 31°07/05.29" 118°00'05.44" 22 6.812 0.812 0.728 6.65 —0.115 
LU21 31°14/28.10" 118°17/15.66" 30 8.375 0.913 0.765 7.65 —0.193 
LU22 31°55/43.05" 118°39/39.40” 11 6.625 0.876 0.758 — —0.155 
LU23 32°09'59.47" 118°49'46.83" 21 6.312 0.82 0.715 6.17 —0.146 
LU24 32°14/35.96" 119°24’44.20" 14 6.812 0.843 0.774 — —0.089 
LU25 31°32'19.00" 120°13'23.19" 12 6.125 0.908 0.745 — —0.22 
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Fig. 2 Geographic patterns of genetic diversity across the species’ range (a) Ho and (b) A,. Piecewise linear correlation identified at 
the boundary point (115.7°E) of the upstream (@) and downstream (A) of Yantgze River. 


Table 2 The analysis of molecular variance (AMova) of Miscanthus lutarioriparius across the species range 


Source df SS MS 

Among Pops 24 941.964 39.248 
Within Pops 616 6177.400 9.980 
Total 641 7119.364 49.228 


and by Dest (above the diagonal). A total of 283 pairwise 
comparisons were significant (P < 0.05) for genetic dif- 
ferentiation between populations after correction for 
multiple tests. LU3, LU8, LU11, LU14, LU18, LU19, 
LU20, LU23 were significantly differentiated from all 
populations after correction for multiple tests (P < 0.01). 
The nonsignificant differentiation was primarily found 
for comparisons within LU17 and LU22 (Table $2), 
which may be due to the limited sample sizes for the 
two populations. Levels of pairwise population differen- 
tiation were low and moderate, ranging between 0.0037 
and 0.1122 (mean Fsr + SEM: 0.063 + 0.008). Consistent 
with these results, values of Jost’s Dest were similar (0- 
0.1993) between all populations (mean Dest + SEM: 
0.153 + 0.02). The aMova showed that the overall popu- 
lation differentiation was low, but significant 
(Ope = 0.102, P = 0.01). Of the total genetic variation 
partitioned in the 25 M. lutarioriparius populations, 10% 
was attributed to the differences among populations 
(df = 24, P < 0.0001), whereas 90% was attributed to the 
differences among individuals (df = 619, P < 0.0001) 
(Table 2). 

STRUCTURE software with the Bayesian model was used 
to estimate the number of homogeneous gene pools (K) 
and confirmed a constant increase in In P(K) from K = 1 
to K = 15 with low variance across replicate runs. Eval- 
uation of the most appropriate K value following the 
procedure by Evanno et al. (2005) found two clear peaks 
for AK, at K = 2 and 4, suggesting that a model with 


Est. Var. Var.% Dye P 
1.139 10% 
9.980 90% 0.102 0.010 
11.119 


two gene pools captures a major split, with substantially 
additional resolution provided under a model with 
K = 4 (Fig. 3a). Considering the inspection of Q (propor- 
tional membership) values for K = 2 with high levels of 
admixture, the K = 4 model was used as prior informa- 
tion. The four clusters of individuals clearly separated 
these collecting locations into stable gene pools. Some 
individuals from downstream locations (LU17-25) had 
more similar alleles than others from the upstream 
populations, with fewer cases in the opposite direction 
(Fig. 3b). Further investigation using hierarchical STRUC- 
TURE analyses revealed that LU8 was split into two 
genetically distinguishable groups, one that had similar 
alleles with other populations and one that was specific 
for this population, suggesting the presence of a barrier 
that partially cut off LU8 from the other populations. 


Spatial genetic structure 


There was a significant association between genetic dif- 
ferentiation and geographic distance across the total 
range, although only a small proportion of the genetic 
variance was explained by geography (Fsr: r = 0.3907, 
P = 0.0002; Dest: r = 0.2947, P = 0.004; Fig. 4). The char- 
acterization of spatial genetic structure at the individual 
level, as calculated with GenAlex, revealed a consistent 
pattern with IBD, where the pairwise kinship coeffi- 
cients decreased with distance. Positive kinship coeffi- 
cients (r) were found in the range of 60 m, with an x- 
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Fig. 3 Estimations of the number of gene pools across the M. lutarioriparius distribution range. (a) The number of Bayesian clusters 
(K) among 25 populations using Ln P(K) (@; Pritchard et al., 2000) and AK (W; Evanno et al., 2005); (b) population structure inferred 
from admixture analysis from K = 2 to K = 7 clusters identified in Bayesian analyses. A thin vertical line represents individuals that 
were separated into K-colored segments. Black lines separate individuals from different populations (labeled above). 


intercept of 78 m (Fig. 5a). The kinship coefficients were 
highest at the first distance class and rapidly decreased 
from there. At the range of 20 m, r value was the high- 
est in the LU2 (r = 0.325). Similarly, autocorrelation in 
all sites was positive up to about 60 m, and, with the 
exception of LU3, LU18, LU23, LU24, and LU25, with r 
values below zero (data not shown). In Fig. 5b, the cor- 
relation values were positive and significant up to 
150 m and significantly negative at 300 m, with an x- 
intercept of 190.8 m. Figure 5a showed an oscillation of 
high and low autocorrelation, indicating a mosaic of 
higher and lower densities at different scales. 


Barrier analysis 


Analyses using Barrier 2.2 identified six main genetic 
discontinuities based on pairwise Fsr values (Fig. 6a). 
The first barrier, computed on the Fsr distance matrix, 
divides LU23 from adjacent populations (LU22 and 
LU24) and was supported completely by 15 loci. The 
results showed that the thickness of the barrier 


decreased from the LU22 to LU24, thus suggesting that 
the genetic barrier was strong only between LU22 and 
LU23. The second break and third break separated 
LU21 and LU22, and LU19 and LU20 with a support of 
11 loci. The fourth break separated samples of LU16 
from LU15 and LU17, which was wholly supported by 
10 loci. This break was in agreement with other analy- 
sis, such as STRUCTURE, which splits these populations 
into upstream and downstream. A fifth break and sixth 
break were detected by 9 loci, of which the former sepa- 
rated LU8 from other neighboring samples and the 
latter separated LU24 from LU25. Other genetic discon- 
tinuities identified by barrier 2.2 did not meet the major- 
ity rule criterion of being supported by at least 50% of 
all loci and therefore was not considered further. Over- 
all, the six geographical barriers significantly and stably 
separated the downstream populations into distinct 
gene pools. When considering genetic structure based 
on the PCoA performed on allelic frequencies, a similar 
pattern was found. PCoA explained 77.73% of the total 
gene variation, on the basis of the first dimension 
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Fig. 4 Mantel test of correlation of genetic distance with geographic distance for populations of M. lutarioriparius. Pairwise values 
between populations are calculated by Fsr (@;Weir & Cockerham, 1984) and Dest (O; Jost, 2008). 
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Fig. 5 Correlogram for the spatial pattern of M. lutarioriparius based on multiple populations in the program GenAlex 6.0. (a) 540 m 
(b) 1000 m. r indicates correlation coefficient. U and L indicate upper and lower limits by 95% confidence interval of the null hypothe- 


sis of no spatial structure. 


(34.05%, P < 0.001), the second dimension (24.07%, 
P < 0.001), and the third dimension (19.61%, P < 0.001) 
(Fig. 6b). 


Gene-flow and population history 


An assignment test showed that individuals from their 
population of origin were 54.3% or higher in all the nat- 
ural populations, and 4.5% or less of the individuals 
were assigned to their nearest origin population. Addi- 
tionally, 37.3% or less of the individuals was assigned 
to the other populations. Some 3.9% or less of the indi- 
viduals were not assigned to any of the populations 
sampled (Fig. S1). Thirty first-generation immigrants 
were detected, of which sixteen had migrated up from 
the downstream populations and fourteen had migrated 
down from upstream (Table 3). Based on Fsr- and 


private allele methods, historical gene flow was 
estimated at 3.718 and 2.862 among M. lutarioriparius 
populations, respectively (Table S1). 

We tested for signatures of recent bottlenecks for 
these populations and found significant heterozygosity 
excess in all populations, except LU6 under the IAM. 
Interestingly enough, none of the populations showed 
significant heterozygosity excess under TPM and SMM. 
The allele frequency distribution tests did not show any 
shifts in allele frequency distribution shapes (Table 4). 


Discussion 


Genetic diversity and adaptation 


For a species that is endemic in wetlands of the Yangtze 
River, M. lutarioriparius exhibited a higher level of 
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Fig. 6 Genetic landscape of M. lutarioriparius. (a) Major discontinuities in gene flow suggested by BARRIER using consensus barriers 
(>8 loci), with rank order given in letters. The A was supported by 15 loci, the B barrier by 11 loci, the C by 10 loci, and the D by 9 
loci. Additional barriers were identified, but were weakly supported (<8 loci) and were largely restricted to the downstream of the 
species’ range. (b) Nonmetric multiddimensional scaling ordination among upsteam (@) and downstream (A) based on Nei’s genetic 


distance (Da; Nei et al., 1983). 


genetic variation across its distribution range (Hg = 
0.682-0.786; A, = 4.74-8.06; Tables 1 and S1) compared 
to Myricaria laxiflora and Zizania latifolia within the same 
habitat and distributing range (Liu et al., 2006; Chen 
et al., 2012). However, the genetic diversity was compa- 
rable with the riparian pioneer tree species Populus nigra 
along the Drome River (France), which had similar dis- 
persal mechanisms (Imbert & Lefèvre, 2003). Generally, 
pioneer species tend to be less subdivided and more 
genetically diverse than species that establish more 
slowly (Vrijenhoek, 2010). Therefore, combined with the 
growing characters and genetic diversity of the species, 
we inferred that M. lutarioriparius might have estab- 
lished quickly and recently. Notably, we found that the 


clonal reproduction of the species was not as wide- 
spread as commonly assumed. We found only three clo- 
nal genotypes in 644 samples, with rhizomes spreading 
no more than 25 m. Above all, M. [utarioriparius showed 
a considerable level of genetic diversity, which may be 
explained by the specific traits of the species, including 
perenniality, self-incompatibility (Fig < 0), weak clonal 
reproduction, and various modes of distribution. 
Genetic diversity plays an important role in the adap- 
tation capabilities of species and their susceptibility to 
environmental changes (Raymond et al., 2013). Measure- 
ments of the relationship between genetic variation and 
adaptive traits have already been conducted in a num- 
ber of species, which showed that species with high 


© 2015 The Authors. Global Change Biology Bioenergy published by John Wiley & Sons Ltd., GCB Bioenergy, doi: 10.1111/gcbb.12278 


10 J. YAN et al. 


Table 3 First-generation immigrants and the direction of 
migration among 25 populations across the distributional 
range; populations divided as upstream and downstream at 
longitude 115.7°E 


Table 4 Wilcoxon tests for recent bottlenecks in the program 
BOTTLENECK. IAM (infinite allele model), TPM (two-phase 
mutational model), SMM (stepwise mutation model); all with a 
normal L-shaped distribution 


Individual Migrated from Direction of migration P value IAM TPM SMM 
LU1 no.10 LU11 Downstream 0.01 LUI 0.00134 0.7057 0.1297 
LU1 no.20 LU13 Downstream 0.001 LU2 0.00629 0.56189 0.1753 
LU2 no.5 LU6 Downstream 0.008 LU3 0.0001 0.97995 0.3225 
LU2no.14 LU6 Downstream 0.001 LU4 0.00008 0.2114 0.93988 
LU3 no.4 LU13 Downstream 0.001 LU5 0.00011 0.4637 0.0832 
LU5 no.1 LU19 Downstream 0.003 LU6 0.07391 0.2312 0.0576 
LU5 no.4 LU24 Downstream 0.006 LU7 0.00011 0.78195 0.19281 
LU6 no.19 LU15 Downstream 0.009 LU8 0.00002 0.00919 0.05768 
LU7 no.14 LU10 Downstream 0.003 LU9 0.00002 0.82089 0.46375 
LU8 no.23 LU9 Downstream 0 LU10 0.00025 0.93988 0.37546 
LU8 no.24 LU9 Downstream 0.003 LU11 0.00029 0.70572 0.52817 
LU9 no.21 LU5 Upstream 0.005 LU12 0.00008 0.15806 0.66855 
LU10 no.10 LU7 Upstream 0.007 LU13 0.00001 0.09344 1 

LU10 no.11 LU3 Upstream 0.009 LU14 0.00021 0.97995 0.27444 
LU11 no.23 LU1 Upstream 0.007 LU15 0.00029 1 0.25223 
LU12 no.24 LU9 Upstream 0.002 LU16 0.00003 0.01825 0.11667 
LU14 no.9 LU15 Downstream 0.008 LU17 0.00029 0.46375 0.82089 
LU15 no.9 LU20 Downstream 0.007 LU18 0.00134 0.11667 0.01309 
LU15 no.22 LU22 Downstream 0.003 LU19 0.00003 0.97995 0.27444 
LU16 no.21 LU7 Upstream 0.004 LU20 0.0005 0.19281 0.66855 
LU17 no.1 LU22 Downstream 0 LU21 0.00004 0.78195 0.66855 
LU19 no.16 LU1 Upstream 0.007 LU22 0.01099 0.56189 0.17535 
LU19 no.21 LU4 Upstream 0.009 LU23 0.00029 1 0.46375 
LU20 no.3 LU4 Upstream 0.002 LU24 0.00107 0.56189 0.97995 
LU21 no.30 LU19 Downstream 0.001 LU25 0.0005 0.25223 0.4332 
LU22 no.2 LU17 Upstream 0.002 

LU22 no.4 LU20 Upstream 0.005 

LU23 no.2 LU5 Upstream 0 a. . A . 

LU23 no.3 LU6 Upstream 0.005 variation was not associated with the longitude of the 
LU25 no.3 LU15 Upstream 0.004 collection sites. Furthermore, genetic variation was not 


level of genetic diversity and gene flow have higher 
phenotypic plasticity (Reed & Frankham, 2003; Ghalam- 
bor et al., 2007; Gruber et al., 2013). For example, LU8 
with lower genetic variation had a weaker adaptability 
than other populations in the Loess Plateau, and other 
populations with higher genetic variation had a stron- 
ger adaptability to establish in a new environment (Yan 
et al., 2012). M. lutarioriparius shows a high level of 
genetic diversity, which is consistent with its capability 
to adapt to the stressful environments. 


Geographic pattern of genetic variation 


It is generally assumed that the effective population size 
and genetic variation increases in downstream popula- 
tions of riparian species due to the movement of seeds 
and propagules with the water. However, M. lutarioripa- 
rius does not seem to follow this pattern. Genetic 


different between the upstream and downstream envi- 
ronment of the Yangtze River when taking the 
traditional longitude breakpoint (115.7°E) based on 
hydrological characteristics. Generally, clinal variation 
in allelic frequencies results from founder events or spa- 
tially restricted gene flow (Vasemägi, 2006; Hasselman 
et al., 2013). The observed pattern therefore indicates 
that M. lutarioriparius populations have enough gene 
flow and are in gene-flow-drift equilibrium. Genetic 
variation slightly decreased in downstream populations 
compared with upstream, which may be due to differ- 
ences in landscape architecture, anthropogenic factors 
such as farming and general habitat suitability. 


Genetic differentiation and gene flow 


Life style, mating system and mode of seed dispersal 
are important adaptive traits that affect genetic differen- 
tiation among populations (Mable & Adam, 2007; Vähä 
et al., 2007; Yan et al., 2009). For instance, perennial, out- 
crossing, bird-dispersed, and wind-disseminated seeds 
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and pollen generally have less genetic differentiation 
than annual, self-crossing species with gravity and 
insect-mediated dispersal. M. I[utarioriparius exhibited 
low genetic differentiation among populations across 
the geographic range in this study, as supported by the 
AMOVA results with only 10% variation existing among 
populations (P = 0.01). Assignment analysis revealed 
that close to 50% of the sampled individuals might 
derive from other populations. Similarly, the analysis of 
gene-flow history also indicated a high rate of migra- 
tion. Gene flow for M. lutarioriparius was estimated to 
be higher than in other terrestrial plants such as Vo- 
uacapoua americana and Castanopsis sclerophylla (Dutech 
et al., 2004; Wang et al., 2011). The low genetic differen- 
tiation and high level of gene flow are likely the result 
of a combination of factors, including a perennial life 
style, self-incompatibility, strong dispersal, and anthro- 
pogenic displacement. 

It is interesting to note that immigration was higher 
in upstream populations than in downstream popula- 
tions. The result may be due to several factors; firstly, 
the interseasonal flood in upstream region may have 
enabled propagules (fruits, seeds, or rhizomes) to dis- 
perse by floating and establishing in wetland after 
floods. Secondly, Dongting lake, as a typical estuary 
wetland, is one of the important stopovers of migration 
route for eastern birds (Deng & Wang, 2000), which 
may have increased the dispersal of caryopses. Finally, 
commercial activity with the plant was higher in the 
upstream regions, which may have increased anthropo- 
genic gene flow. 


Spatial genetic structure, population history, and 
population connectivity 


The detection of four genetically distinguishable 
groups among the 25 populations by Bayesian cluster- 
ing should be considered conservative, as this method 
has limited power under scenarios involving weak 
genetic differentiation (Faubet et al., 2007; Faubet & 
Gaggiotti, 2008). Despite the low population differenti- 
ation, we found a general pattern of IBD at the scale 
of kilometers, which indicates either ongoing gene 
flow or a recent colonization process (Slatkin, 1993). 
Additionally, kinship coefficients between individuals 
were positive at a scale of hundreds of meters and 
only disrupted over distances of around 200m by 
spatial autocorrelation analysis, a widespread tool for 
estimating genetic structure and patterns of gene flow. 
The correlation between genetic distance and geo- 
graphic distance is consistent with population mixing 
by migratory movements. The significant isolation-by- 
distance patterns and the spatial autocorrelation pat- 
terns are indicative that the populations were close to 


equilibrium. Moreover, the significant isolation of dis- 
tance allowed us to reject the infinite island model of 
migration, suggesting that the species was not sub- 
jected to bottleneck effects. 

Detecting genetic discontinuities is necessary to 
understand dispersal restrictions and to identify factors 
which influence population structure (Manel et al., 2003; 
Meeuwig et al., 2010). Several genetic discontinuities 
were detected across the distribution range of M. luta- 
rioriparius (Fig. 6a), despite the high dispersal capacity 
of the species. PCoAs further supported the pattern 
(Fig. 6b). We found that barriers to individual migra- 
tions largely occurred in regions with high anthropo- 
genic activities (e.g., dam and urban development). The 
three major genetic borders (LU8, LU16, and LU23) cor- 
respond with physical barriers, such as LU8 along the 
villages, and LU16 and LU23 along the Yangtze River 
levees. Similar patterns have been found in other plant 
species. Small irrigation diversion dams, for instance, 
were found to restrict the interactions of Oncorhynchus 
mykiss populations (Weigel et al., 2013). 

In conclusion, this study assessed genetic variation 
and gene flow across M. lutarioriparius populations. Our 
results provide a genetic basis for the adaptation 
potential of the species that was revealed by some pre- 
vious studies (Yan et al., 2012, 2015; Xu et al., 2015). M. 
lutarioriparius has a high level of genetic diversity and a 
low level of genetic differentiation, both of which can 
enhance the capacity of plants to adapt to changing 
environments. Furthermore, the species shows long-dis- 
tance, bidirectional gene flow with unequal dispersal 
among populations in the middle and lower Yangtze 
River, probably as a result of anthropogenic activities 
and landscape features. The high level of gene flow that 
in turn enriched the genetic variation of the species 
should be beneficial for the domestication and cultiva- 
tion of the energy crop on marginal land with unfavor- 
able and often stressful conditions. 
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